Image registration of multiple medical imaging modalities using a multiple degree-of-freedom-encoded fiducial device

ABSTRACT

Disclosed is a system and method for registering images from two medical imaging modalities in a six-degree-of-freedom coordinate space, for the purpose of providing real time registered imagery for optimization of medical procedures. The system and method uses fiducial that is visible to the first imager and in a fixed position and orientation relative to the second imager. The first imager may be an X-ray device such as a C-arm fluoroscope, and the second imager may be a Transrectal Ultrasound (TRUS) device.

This application is a divisional of U.S. application Ser. No. 11/998,452filed Nov. 30, 2007, which is a continuation of U.S. application Ser.No. 10/895,398 filed Jul. 21, 2004, which claims priority to U.S.Provisional Application No. 60/488,965 filed Jul. 21, 2003, the entirecontents of these applications are incorporated herein by reference.

Research and development efforts associated with the subject matter ofthis patent application was supported by the National Science Foundationunder Grant No. ERC 9731478, and by the National Institutes of Healthunder Grant No. 1R43CA099374-01.

BACKGROUND OF THE INVENTION

1. Field of the Invention

The present invention relates to the registration of multiple medicalimaging modalities in support of real time optimization of medicalprocedures. In particular, the invention relates to the registration ofultrasound and C-arm fluoroscopy imagery for these purposes.

Adenocarcinoma of the prostate is the most commonly diagnosed cancer inthe U.S. male population. During the last half decade, there have beenapproximately 200,000 new cases of prostate cancer diagnosed each year,which is comparable to breast cancer diagnosis, and there is no evidencethat this number would significantly decrease in the foreseeable future.For several decades, the definitive treatment of low-risk prostatecancer was radical prostatectomy or external beam radiation therapy(EBRT). During the 90's the technique of transrectal ultrasound (TRUS)guided transperineal low dose-rate brachytherapy became a well-proventreatment alternative, comparable to surgery and EBRT, has demonstratedexcellent local control rates, and most recent data reports excellentlong-term (10-12 year) disease—free survival rates equivalent to radicalprostatectomy and external beam radiotherapy. Moreover, studies indicatethat there are no unexpected urethral or rectal complications or sideeffects with low dose-rate brachytherapy.

However, the rate of rectal and urethral complications associated withlow dose-rate brachytherapy is still high. The underlying reason forthis has predominantly been the lack of adequate visualization andcontrol of the implant process, which leads to improper placement of theimplanted radiation source and inaccurate delivery of the sources ofradiation dosage. As one skilled in the art will readily appreciate,radioisotopes are used for non-invasive visualization and controlprocesses such as low dose-rate brachytherapy. There have been advancesin the ability to use different isotopes to achieve better dosedistributions. However, even the best isotope is useless withoutadequate visualization, localization, and subsequent control of theimplanted radiation source. The ability to intraoperatively localizeimplanted radiation sources (hereinafter “seeds”) relative to theprostate is key to enabling dynamic dose calculation during theprocedure. The solution of this visualization and control problem wouldreduce the probability of inaccurately placed implants, thus presentingan opportunity for improved outcomes.

TRUS imaging generally provides satisfactory differentiation of relevantsoft issue, but implanted brachytherapy seeds cannot be clearlyidentified in the TRUS images. Advancements in ultrasound equipmenttechnology are expected to facilitate identification and localization ofthe seeds in the future, but such equipment will not be available forpractitioners in the foreseeable future. On the other hand, currentlysixty percent or more of the practitioners use intra-operative C-armfluoroscopy for purposes of visualization, localization, and control.Although C-arm fluoroscopy, which employs transluminal X-Ray imagery,can accurately localize seeds, it does not differentiate soft tissue.Hence, there is a strong demand for using C-arm fluoroscopy togetherwith a TRUS-guided delivery system.

2. Discussion of the Related Art

Problems and limitations inherent in the related art may be properlyillustrated by discussing a particular combination of imaging modalitiesusing C-arm X-ray fluoroscopy and transrectal ultrasound (TRUS). C-armX-ray fluoroscopy may be the most widely used intra-operative imagingmodality in general surgery, and approximately 60+% of the prostatebrachytherapy practitioners use it for qualitative implant analysis inthe operating room in non-computational qualitative manner. While C-armfluoroscopy has been used intra-operatively, it has not been utilized inquantitative intra-operative analysis. C-arm has been used as a sologuidance modality. However, since TRUS emerged as a primary imageguidance modality, C-arm fluoroscopic x-ray imaging has become asecondary tool for gross visual observation. Very few attempts have beenmade to relate fluoroscopic images to soft tissue anatomy, and withlittle success. These attempts have generally used thin metal wireinside a Foley catheter to visualize the prostatic urethrafluoroscopically in anterior-posterior and lateral projections. In otherapproaches, gold marker seeds have been implanted into the prostate, andthe relative positions of the needles and marker seeds have beenobserved in fluoroscopy.

X-ray radiography has been used extensively for post-implantbrachytherapy evaluation using multi-view X-ray to recover seedlocations post implant and determine gross dosimetry. Here thefundamental problem is matching a large number of seeds with theirprojections in multiple X-ray images when some seeds obscure each otherand solid objects also can get in the way. Automated methods have beenexplored, but they also assume conditions that most likely cannot be metwith a C-arm fluoroscope, particularly in a realistic intra-operativescenario. Such assumptions include no extrinsic object in the field,optimal beam energy, arbitrary number and orientation of X-ray shots, orunlimited processing time and computational resources—none of which isrealistic in real-time image-guided surgery. In conclusion,“off-the-shelf” post-implant seed matching and reconstruction techniquesof the related art cannot be expected to work in the operating room.Many approaches have met with great difficulty due to the inability toaccurately determine the imaging angles relative to the prostate forreconstruction of the multiple projection images. This not only reducesaccuracy, but makes the process too lengthy to be used intraoperatively.

The use of implanted needles as fiducial markers for registration ofbiplane TRUS data has been explored, but several key problems have beenleft unsolved. First, the use of implanted needles as fiducials may notbe practical, because most practitioners implant only one needle at atime and they do not use stabilizing needles. Second nearly paralleltransperineal needles, when used as fiducials, encode weakly in theapex-base direction, with little spatial resolution. Third, becauseX-ray and TRUS imaging are not simultaneous, it is imperative that thefiducials do not move relative to the prostate until both imagingsessions are complete. For example, during TRUS scanning, the prostatedeforms, dislocates, and the needles dislocate relative to the prostate.Fourth, previous attempts at using implanted needles as fiducial markersdid not account for the need to pre-operatively calibrate the C-armfluoroscope, including removing image distortion and some form ofintra-operative tracking to know where the multiple X-ray shots arecoming with respect to one another.

Accordingly, TRUS guided transperineal low dose-rate brachytherapy hasemerged as one of the preferred treatments for low-risk prostate cancer.While ultrasound is generally an excellent tool for guiding the implantneedles with respect to prostate anatomy, it cannot reliably show thelocation of radioactive seeds after they are released in the prostate.Intraoperative C-arm fluoroscopy does adequately show the location ofimplanted seeds, but it cannot differentiate soft tissue; thus, itcannot be used for visualizing prostate anatomy.

Accordingly, what is needed is intra-operative fusion of these twocomplementary modalities, which offers significant clinical benefit byallowing for real-time optimization of the brachytherapy implant as theprocedure progresses in the operating room.

The problems discussed are not limited to TRUS guided transperineal lowdose-rate brachytherapy discussed above, but apply to proceduresrequiring accurate placement of objects such as surgical tools orimplants within a patient. The problems and limitations apply to anyprocedure in which two imaging modalities are used, and in which thereis a need to fuse or merge the products of the two imaging modalitiesinto a single image space.

SUMMARY OF THE INVENTION

Accordingly, the present invention is directed to the image registrationof multiple medical imaging modalities using a multipledegree-of-freedom encoded fiducial device that substantially obviatesone or more of the problems due to limitations and disadvantages of therelated art.

An advantage of the present invention is to provide improvedvisualization and control during invasive surgical procedures.

Another advantage of the present invention is improved localization ofsurgically implanted objects in real time.

Another advantage of the present invention is to enable real timeregistration of multiple medical imaging modalities, whereby eachmodality provides images of different objects of distinct compositionsin the same volume.

Another advantage of the present invention is to better enable dynamicdose calculation during radiation implant procedures, including theability to adjust the planned implant locations to better optimizedosage as the implant progresses.

Additional features and advantages of the invention will be set forth inthe description which follows, and in part will be apparent from thedescription, or may be learned by practice of the invention. Theobjectives and other advantages of the invention will be realized andattained by the structure particularly pointed out in the writtendescription and claims hereof as well as the appended drawings. Toachieve these and other advantages and in accordance with the purpose ofthe present invention, as embodied and described, a system forregistering images from two medical imaging modalities comprises a firstimager having a first coordinate frame; a second imager having a secondcoordinate frame; a fiducial having features visible to the firstimager, the fiducial having a fiducial coordinate frame that has a knownposition and orientation relative to the second coordinate frame; and adata system.

In another aspect of the present invention, the aforementioned and otheradvantages are achieved by a device for registering images acquired by aplurality of medical imaging modalities. The device comprises a supportstructure; a first wire disposed on the support structure, the firstwire substantially formed in the shape of a straight line; and a secondwire disposed on the support structure, the second wire substantiallyformed in the shape of an ellipse.

In another aspect of the present invention, the aforementioned and otheradvantages are achieved by a method for registering images from a firstimaging modality to a second imaging modality. The method comprisesacquiring a first image from a first imager having a first coordinateframe; acquiring a second image from a second imager having a secondcoordinate frame; determining a first coordinate transformation betweenthe first coordinate frame and a fiducial coordinate frame correspondingto a fiducial; and registering the first image to the second coordinateframe according to the first coordinate transformation

In another aspect of the present invention, the aforementioned and otheradvantages are achieved by a method for reconstructing athree-dimensional structure of a plurality of objects using imagesacquired by a two-dimensional imager. The method comprises acquiring afirst image of the plurality of objects at a first orientation;estimating a first pose corresponding to the first orientation;acquiring a second image of the plurality of objects at a secondorientation; estimating a second pose corresponding to the secondorientation; acquiring a third image of the plurality of objects at athird orientation; estimating a third pose corresponding to the thirdorientation; computing a first minimum cost correspondence between eachof the seeds in the first image with each of the seeds in the secondimage and with each of the seeds in the third image; and computing afirst minimum cost pose estimation based to the first correspondence.

It is to be understood that both the foregoing general description andthe following detailed description are exemplary and explanatory and areintended to provide further explanation of the invention as claimed.

BRIEF DESCRIPTION OF THE DRAWINGS

The accompanying drawings, which are included to provide a furtherunderstanding of the invention and are incorporated in and constitute apart of this specification, illustrate embodiments of the invention andtogether with the description serve to explain the principles of theinvention.

In the drawings:

FIG. 1 illustrates a first exemplary system for registering imagery froma first imaging modality to that of a second imaging modality accordingto the present invention;

FIG. 2 shows two distinct images of a fiducial, as may be acquired by aC-arm fluoroscope at multiple angles;

FIG. 3 illustrates an exemplary process for registering imagery from aC-arm fluoroscope to that of a TRUS imager according to the presentinvention;

FIG. 4 illustrates an exemplary “basic” process for estimating the poseof a C-arm fluoroscope relative to a fiducial;

FIG. 5 illustrates an exemplary “robust” process for estimating the poseof a C-arm fluoroscope relative to a fiducial;

FIGS. 6A and 6B show a field of brachytherapy seeds as acquired by aC-arm fluoroscope at multiple angles;

FIG. 7 illustrates multiple C-arm fluoroscope fields of view capturing asingle field of brachytherapy seeds from multiple viewing angles;

FIG. 8 illustrates an exemplary process for matching brachytherapy seedsin multiple registered images according to the present invention;

FIG. 9 illustrates a second exemplary system for registering C-armfluoroscope images to TRUS images according to the present invention;

FIG. 10 illustrates an exemplary process for registering ultrasound andfluoroscopy according to the present invention;

FIG. 11 is another illustration of an exemplary process for registeringultrasound and fluoroscopy according to the present invention;

FIG. 12 a shows an exemplary fluoroscope dewarping kit as may be used inthe present invention;

FIG. 12 b shows an exemplary dewarping kit installed on a C-armfluoroscope;

FIG. 13 illustrates an exemplary transrectal sheath, along with anultrasound probe, as may be used in the present invention;

FIG. 14 shows an exemplary transrectal sheath mounted to a TRUS stepperbase;

FIGS. 15A-15C are cutaway views of a transrectal sheath, a TRUS probe,acoustically coupled according to the present invention;

FIG. 16 shows a prostate phantom that may be used in a pre-opcalibration step of the present invention;

FIGS. 17A-17B show an exemplary transrectal sheath, which includes airchannels, an X-ray fiducial, and an acoustic fiducial;

FIG. 18 illustrates an exemplary system for performing seedreconstruction according to the present invention;

FIG. 19 illustrates three C-arm images of a plurality of seeds acquiredat different poses; and

FIG. 20 illustrates an exemplary process for registering objects in a3-D image space using images acquired at multiple poses with a single2-D medical imaging system.

DETAILED DESCRIPTION OF THE ILLUSTRATED EMBODIMENTS

The invention provides registration of two medical imaging modalities byregistering multiple images from a first imaging system, preferably atdifferent angles, to a three-dimensional coordinate frame defined by thefield of view of a second medical imaging system. The two medicalimaging modalities preferably complement each other such that oneimaging modality provides better quality images of certain features thanthe other, whereas the second imaging modality provides better qualityimages of other features. Registering the images from the two imagingmodalities is done by determining a pose estimation, wherein the poseestimation is a six degree of freedom (six-DOF) coordinatetransformation. The six-DOF transformation may be represented as atransformation matrix, or some other mathematical representation, fromthe coordinate frame of a fiducial, such as a fiducial device, to thecoordinate frame of the field of view of the first imaging modality. Thefiducial device is a six-DOF spatially encoded reference that is visibleto the first imaging modality and has a known six-DOF coordinatetransformation relative to the coordinate frame of the second imagingmodality. Registration, as used herein, is the act of using the poseestimation and the known six-DOF transformation matrix from thecoordinate frame of the fiducial device to the coordinate frame of thesecond imaging modality to map the image from the first imaging modalityinto the coordinate frame of the second imaging modality. This mappingmay be done of a pixel by pixel basis. By registering the imageryprovided by the two imaging modalities in such a manner, the inventionprovides a user, such as operating room personnel, with the capabilityof locating an object in three-dimensional space within the patient inreal time. Examples of this capability include being able to locateradiation seeds injected into a patient's prostate while knowing thelocation of the seeds relative to soft tissue within and surrounding theprostate.

FIG. 1 illustrates a first exemplary system 100 for registering a firstand second medical imaging modalities in a real time according to thepresent invention. System 100 includes a first medical imaging subsystem115 (hereinafter “first imager”), which has a first imager imageprocessor 130; a second imaging subsystem 110 (hereinafter “secondimager”), which has a second imager image processor 135; a fiducialdevice 105; a data system 140, and a user interface 145.

First imager 115 and second imager 110 may be selected based onwavelength regimes (such as X-ray), electromagnetic properties (such asMRI), or acoustic properties (such as ultrasound), whereby the firstimager 115 may be suited for imaging objects such as metal instrumentsor radioactive seeds, and second imager 110 may be suited for imagingsoft tissue.

The subjects of the first imager 115 and the second imager 110 mayinclude a region of soft tissue 120 and foreign objects 125 introducedinto the soft tissue 120. Examples of soft tissue 120 may include agland such as a prostate, or an organ. Examples of foreign objects 125may include sensors, radioactive seeds, surgical implants, instrumentssuch as needles or catheters, or some object or objects that have acomposition different from that of the soft tissue 120.

First imager 115 may be a C-arm fluoroscope, although other transluminarprojective imaging systems may be used that provide imagery frommultiple angles, provided that the first imager 115 is capable ofproviding images of the foreign objects 125 and the fiducial device 105with sufficient spatial resolution and quality to allow the orientationof the fiducial device 105 and the positions of the foreign objects 125to be determined to a desired accuracy and precision. The first imager115 may be repositioned between image acquisition. Alternatively, if thefirst imager 115 simultaneously provides multiple images from differentfields of view, such as a first imager 115 having multiple spatiallycalibrated detectors, then the first imager 115 may acquire multi-angleimagery without having to be repositioned. Other examples of candidatefirst imagers 115 include X-ray Biplane simulators, X-ray Angiography,MV imaging, infrared telemetry, Nuclear Imaging like PET (PositronEmission Tomography). It will be readily apparent to one skilled in theart that still other medical imaging systems may be used for the firstimager 115.

The second imager 110 provides images having a coordinate system towhich the images from the first imager 115 are registered. The secondimager 110 preferably provides images in a three-dimensional space.Second imager 110 may be a Transrectal Ultrasound (TRUS) system,although other alternate imaging systems may be used that provide athree-dimensional image to which imagery from the first imager can beregistered. Other possible imaging systems for second imager 110 includecomputed tomography (CT), single photon emission computed tomography(SPECT), positron emission tomography (PET), kilovoltage/megavoltage(KV/MV) imaging, optical imaging, ultrasound, and MRI (MagneticResonance Imaging) systems. Preferably, the second imager 110 providesreal time imagery of soft tissue 120. Although the second imager 110 ofexemplary system 100 includes a TRUS system that is inserted into a bodycavity, other systems may be used that have sensors placed outside thepatient's body and which may or may not make contact with the patient'sbody.

The fiducial device 105 is a spatial reference designed to encode sixdegrees of freedom (DOF) information from a single image from any angle.Six-DOF includes three axes of translation and three angles of rotationaround the respective three axes, wherein the three axes define thecoordinate frame of the fiducial device 105. The fiducial device 105 hasfeatures 106 that are visible to the first imager 115. The features 106have a specific shape and orientation, and are made of a material thatis visible to the first imager 115. Encoding six DOF in an image fromany angle means that the features 106, as imaged by the first imager115, yield a unique projection or shape in the image for any givenviewing angle. In other words, there is a one-to-one mapping between theprojected image of the features 106 disposed on the fiducial device 105and the viewing angle of the first imager 115.

FIG. 2 illustrates how the design of the features 106 disposed on thefiducial device 105 may yield a unique image depending on viewing angleof the first imager 115. Images 210 show the features 106, which aredesigned to be visible to the first imager 115. Each mage 210 of thefeatures 106 on the fiducial device 105, as seen by the first imager115, shows a unique projection of the fiducial device 105 as a functionof viewing angle.

The accuracy, precision, and robustness of pose estimation depends onthe design and precise manufacturing of the fiducial device 105.Fiducial device 105 may have a group of features 106 such as ellipses,straight lines, and points disposed on the device. The choice offeatures are not limited to the above. Other exemplary features 106 mayinclude parabolas, hyperbolas, helixes or other parametric curves.

The use of ellipses, straight lines, and points provide differentadvantages and disadvantages. An advantage of using a givenconfiguration of features 106 is that it may provide for an accuratepose estimation over a wide range of orientations. In other words, theshape of the given feature may be such that there is a high probabilitythat resulting pose estimation will be accurate. A disadvantage may bethat the shape is complex and pose estimation is very computationallyintensive. In short, the shape of a given feature within features 106may provide a robust pose estimation in that it provides for a nominallyaccurate pose estimation for a short and consistent amount of computingtime. For instance, ellipses, straight lines, and points respectivelyincrease in order of possibility of error of 6-DOF determination, butrespectively decrease in order of computational complexity. In otherwords, ellipses make pose estimation more accurate, while points makepose estimation more robust. In the embodiment illustrated in FIGS. 1and 3 fiducial device 105 has features 106 that include a combination ofshapes that provide a balance of such advantages and disadvantages,thereby substantially balancing accuracy and robustness.

For any feature within features 106, regardless of its shape, theaccuracy of pose estimation is a function of orientation. For example,the accuracy of pose estimation using a single ellipse is greater from aviewing orientation normal to the plane defined by the ellipse, andlesser from a viewing angle parallel to the plane. Accordingly, byproviding a particular combination of shapes and orientations forfeatures 106 on fiducial device 105, fiducial device 105 may be designedsuch that the accuracy of a given pose estimation is consistent atdifferent viewing angles. For instance, the orientation of features 106should be such that with a change in first imager 115 orientation, asthe orientation of one feature within features 106 in the image becomesless discernable, the orientation of another feature within features 106should become more discernable. With such a configuration of features106, the accuracy and precision of the six-DOF encoding of fiducialdevice 105 may be minimally dependent on its orientation, which meansthat the precision of the pose estimation is consistent regardless ofthe orientation of fiducial device 105. The combination of features 106on the fiducial device 105 illustrated in FIG. 1 is exemplary, and othercombinations are possible and within the scope of the invention.

Fiducial device 105 may be fixed in 6-DOF space relative to the field ofview of the second imager 110, which would substantially fix thecoordinate frame of the fiducial device 105 to the coordinate frame theimage space of second imager 110. Optionally, fiducial device 105 may beable to move relative to the field of view of second imager 110, inwhich case there would need to be a way of measuring the position andorientation (six-DOF) of fiducial device 105 relative to second imager110. Measuring the position and orientation may be accomplished bydevices such as optical encoders, mechanical encoders, and the like. Thesix-DOF measurements made by such encoders may then be used to derive atransformation matrix from the coordinate frame of the fiducial device105 to the coordinate frame of the second imager 110. Fiducial device105 may be integrated with, co-located with, or located separately fromthe second imager 110. Fiducial device 105 may be placed within thepatient, such as within a cavity, or may be placed outside the patient,provided that fiducial device 105 remains substantially fixed relativeto the tissue 120 and objects 125 within the tissue.

Features 106 of fiducial device 105 may be formed out of any combinationof substances that are appropriately opaque to first imager 115. Forexample, fiducial device 105 of exemplary system 100 includes features106 formed of about 1 mm thick Tantalum wire. The Tantalum wire isdisposed on a hollow cylinder, which serves as a support structure.Fiducial device 105 of exemplary system 100 may have dimensions suchthat the features 106 fit within a volume of about 3 cm×3 cm×5 cm. Iffirst imager 115 is X-ray based, then features 106 are preferably formedof materials with good visibility under X-ray illumination, such asTantalum. Other materials may be used in the case of an X-ray-basedfirst imager 115, provided that the material results in a good X-raycontrast. If first imager 115 operates in a different imaging modality,then features 106 should be made of material that will providesufficient contrast in the images of first imager 115.

Data system 140 includes at least one processor 141 for executingcomputer instructions associated with the present invention; a memory142 encoded with computer instructions and configuration constants; andinterface devices 143 for communicating with the first imager imageprocessor 130, the second imager image processor 135, and the userinterface 145. Processor 141 may include a single computer or may beinclude multiple computers that may be collocated or distributed over awide area and networked together. Data system 140 may be integrated intoeither first imager image processor 130, second imager image processor135, or distributed between both. If data system 140 is distributedbetween both image processors, then first imager image processor 130 andsecond imager image processor 135 may communicate directly with eachother and the user interface 145. The user interface 145 may becollocated with the data system 140, integrated with the data system140, or may be remotely located. If remotely located, user interface 145may communicate with the data system 140 over the internet. It will bereadily apparent to one of ordinary skill that many such architecturesfor data system 140 are possible.

As mentioned above, the data system 140 includes a memory 142 encodedwith computer instructions and data values. The computer instructions(hereinafter the “software”), when executed by the processor 141,perform the steps to register multiple medical imaging modalities.

FIG. 3 illustrates an exemplary process 300, which may be implementedall or in part by the software, for performing image registration andlocating a plurality of foreign objects 125 within soft tissue 120.Exemplary process 300 may be fully automated, or it may involve operatorinteraction through the user interface 145. Exemplary process 300 may beused with an embodiment of exemplary system 100 in which the firstimager 115 is a C-arm fluoroscope; the second imager 110 is a TRUSsystem; and the foreign objects 125 are brachytherapy seeds that havebeen injected into prostate soft tissue 120.

Exemplary process 300, as implemented using exemplary system 100, worksas follows. In step 305, first imager 115 (the C-arm fluoroscope) ispositioned to a first orientation relative to the fiducial device 105.This may be done manually or automatically through the software executedby the processor 141.

In step 310, first imager image processor 130 acquires fluoroscopy datafrom the first imager 115 and processes the data corresponding to animage. In the case in which first imager 115 is a fluoroscope and secondimager 110 is a TRUS probe, the TRUS probe may be removed from the fieldof view of the first imager 115 to prevent interference with the imageacquired by first imager 115. Data system 140 retrieves the C-arm imagedata from first imager image processor 130 and stores the C-arm imagedata values (hereinafter the “acquired C-arm image”) in memory 142.Further, second image processor 135 acquires TRUS data from the secondimager 110. Data system 140 retrieves the TRUS image data from secondimage processor 135 and stores the TRUS image data values (hereinafterthe “acquired TRUS image”) in memory 142.

Once data system 140 has stored the image data values, processor 141executes the software to compute the six-DOF pose estimation of thefirst imager 115 in step 315. As mentioned earlier, the pose estimationis a 6-DOF representation, such as a transformation matrix, of thetransformation from the coordinate frame of fiducial device 105 to thecoordinate frame of first imager 115.

In step 320, data system 140 executes the software to register theacquired C-arm image (from first imager 115) to the coordinate frame ofthe acquired TRUS image (from second imager 110), by using the poseestimation computed in step 315, and stores the resulting data valuescorresponding to the registered C-arm image (hereinafter the “registeredC-arm image”) in memory 142. Algorithms for performing the registrationof step 320 are well known to the art of image processing. It will beapparent to one skilled in the art that existing image processingsoftware packages may be used to implement image registration in step320, and that the use of such software packages is within the scope ofthe invention.

In step 325, data system 140 retrieves the registered C-arm image datavalues, and the acquired TRUS images, from memory 152, and sends the twoimages to user interface 145 to be displayed in an overlaid fashion.

In step 330, data system 140 executes the software to query the user ifanother registered C-arm image, from first imager 115, is desired. Atleast two registered C-arm images are required to perform seed matchingin subsequent step 340. If the answer is affirmative, steps 305-325 maybe repeated. During a second iteration of steps 305-325, a second TRUSimage need not be acquired by second imager 110, if multiple C-armimages are to be registered to the coordinate frame of a single TRUSimage.

In step 335, data system 140 retrieves the acquired TRUS image and theregistered C-arm images from memory 142.

In step 340, data system 140 executes the software to perform seedmatching, in which the seeds imaged in one registered C-arm image arematched to the corresponding seeds in other registered C-arm images. Inother words, in seed matching, a given seed is identified in bothregistered C-arm images, and then located in the space defined by thecoordinate frame of the second imager 110.

Exemplary processes for implementing step 315 are described as follows.

FIGS. 4 and 5 illustrate two exemplary implementations of step 315. FIG.4 illustrates a “basic” pose estimation algorithm 400 that may beimplemented in the software in data system 140; and FIG. 5 illustrates a“robust” pose estimation algorithm 500 that may be implemented in thesoftware. Selection of which pose estimation algorithm to use may dependon the available computational power of data system 140 as well as poseestimation accuracy and precision requirements.

Referring to FIG. 4, “basic” algorithm 400 is an iterative process inwhich a pose estimation is iteratively computed until it converges on asolution that is within a predefined acceptable error threshold.

In step 410 of “basic” algorithm, processor 141 executes the software,which segments features 106 of fiducial device 105. To segment features106 is to distinguish features 106 of fiducial device 105 from thebackground image. The software may implement any of a number ofsegmentation image processing algorithms that are capable of identifyingthose pixels that correspond to features 106 and identifying parametriccurves corresponding to features 106 to sub-pixel precision. Typicalsegmentation methods may include thresholding followed by edge filteringand Hough Transforms. Approaches like Random Sample Consensus (RANSAC)can be used to further improve results. Various segmentation algorithmssuch as these are known to the art, and their use is within the scope ofthe invention.

The software may implement segmentation in a semi-automatic mode or anautomatic mode. In semi-automatic mode, the software displays theacquired C-arm image (which shows features 106) on user interface 145.The operator then selects various points in the acquired C-arm image(for example, by clicking on it with a mouse attached to user interface145) in a pre-set order. The selected points should be close to theshown features of interest, which are the displayed features 106 in theacquired C-arm image. The software does the following: it performs aminimum intensity centroid search in the neighborhood of each selectedpoint, which fine-tunes the location; it performs a least square curvefitting in order to obtain the equation of the straight line andellipse; and it projects the centroid on its corresponding straight lineand ellipse to get the final position.

In automatic mode, the software implements a Hough Transform to segmentthe parametric curves of the corresponding to features 106 on fiducialdevice 105. The software achieves this by creating an edge image basedon the image acquired by first imager 115, and then searching in theparametric space defined by the edge image for parametric curves, whichincludes straight lines and ellipses, that correspond to features 106 offiducial device 105.

To identify straight lines, the software performs a Hough Transform. Thesoftware may then fine tune each straight line by doing the following:searching for minimum intensity points in a small neighborhoodsurrounding each line; running a RANSAC-based least square fitting onthese points; and selecting the three strongest near-parallel straightlines. RANSAC-based line fitting includes running a filter in theneighborhood of each line and extracting a large number of possiblepoints on the line. Some of these points may be erroneous and maycorrupt the final answer, while most of the points may be good. Thisprocess may be repeated many times to pick points from the output of theHough Transform at random, with which it fits a line. The line thatrepeats most often is the final correct answer.

To identify ellipses, the software may execute a Hough Transform toestimate various possible axes of symmetry in the image. To achievethis, the software performs a 2D Hough Transform to find the bestellipse for each center of axis of symmetry. The software thenfine-tunes the ellipse associated with each center using an approachsimilar to the one used for the lines. Finally, the software selects thetwo strongest ellipses, which correspond to the ellipses in features 106on fiducial device 105.

The result of the segment fiducial step 410 is a segmented projection offeatures 106 in the two-dimensional coordinate frame of first imager115. The segmented projection may include a set of points in thetwo-dimensional coordinate frame of first imager 115, a set of lineseach defined by two endpoints and parameters defining the line betweenthe endpoints, a set of ellipses defined by a plurality of pointsdefining a given ellipse, and a parametric equation for the ellipsedefined by the points. The segmented projection may be used as aninitial estimate 415 for the pose estimation.

The software then computes the pose estimation by iteratively performingsteps 425-445. In step 425, the software retrieves from memory 142 thedata values corresponding to the current pose estimation and parametervalues corresponding to features 106. As mentioned earlier, if this isthe first pass through step 425, then the current pose estimation isinitial estimate 415. Then the software generates an estimated or “fake”image of features 106, which is a two-dimensional projection of thefeatures 106 according to the current pose estimation. The software thenstores the estimated image data values in memory 142. Alternatively, theprocess of generating an estimated or fake image may be done usingDigitally Reconstructed Radiography.

In step 435, the software retrieves the data values corresponding to theC-arm image acquired by first imager 115, and the estimated or “fake”image, and computes a Euclidean distance or Euclidean error function.The Euclidean distance is a matrix containing vectors between thecorresponding points, lines, and ellipses of the acquired and estimatedimage. The Euclidean distance δx_(i), δy_(i) is a measure of theaccuracy of the pose estimation and is defined by the relation

δx _(i) =x _(i)(R ^(i) , T ^(i))− x _(i)

δy _(i) =y _(i)(R ^(i) ,T ^(i))− y _(i)

where (x_(i), y_(i)) are the coordinates for a given point on theestimated image, ( x _(i), y _(i)are the coordinates of thecorresponding point on the actual image, (R^(i), T^(i)) are therespective rotational and translational components of the current poseestimation, and i corresponds to the ith iteration of steps 425-440. Theabove equations for the Euclidean distance δx_(i), δy_(i) are for pointsin the images. The Euclidean distance for two lines may be computed bycomputing the Euclidean distances for the respective endpoints of eachline. The Euclidean distance for two ellipses may be computed using thefollowing relation:

$\frac{ab}{2\sqrt{a^{2} + b^{2}}}{{{\frac{\left( {{A\; C} - {B^{2}/4}} \right)}{\det (M)}\begin{bmatrix}{\overset{\_}{x}}_{i} & {\overset{\_}{y}}_{i} & 1\end{bmatrix}}\begin{bmatrix}A & {B/2} & {D/2} \\{B/2} & C & {E/2} \\{D/2} & {E/2} & F\end{bmatrix}}\begin{bmatrix}x_{i} \\y_{i} \\1\end{bmatrix}}$

where a is the major axis of the ellipse, b is the minor axis of theellipse, ( x _(i), y _(i)) are the coordinates of a point on the ellipsein the acquired image, (x_(i), y_(i)) are the coordinates of a point onthe ellipse in the estimated image, M is shorthand for the matrix

$\quad\begin{bmatrix}A & {B/2} & {D/2} \\{B/2} & C & {E/2} \\{D/2} & {E/2} & F\end{bmatrix}$

and A-F are parameters defining an ellipse in the acquired image, asdetermined in segmentation step 410. The parameters A-F define anellipse according to the general equation for a conic representation ofan ellipse

${{\begin{bmatrix}x & y & 1\end{bmatrix}\begin{bmatrix}A & {\; {B/2}} & {D/2} \\{B/2} & C & {E/2} \\{D/2} & {E/2} & F\end{bmatrix}}\begin{bmatrix}x \\y \\1\end{bmatrix}} = 0.$

In step 430, the software generates a Jacobian Matrix based on thecurrent pose estimation. In general, the Jacobian Matrix corresponds toa multidimensional derivative of the current pose estimate. In anexemplary implementation, the Jacobian is a 48×6 matrix, where the “48”dimension corresponds to 9 points leading to 18 of the dimensions (x andy), 3 lines leading to 6 dimensions, 2 ellipses leading to 24dimensions, and the “6” dimension refers to the six-DOF.

In an exemplary embodiment, the Jacobian Matrix for a given point may bea combination of translational and rotational elements:

$\begin{matrix}{J = \begin{bmatrix}\frac{\partial x_{i}}{\partial T_{1}} & \frac{\partial x_{i}}{\partial T_{2}} & \frac{\partial x_{i}}{\partial T_{3}} & \frac{\partial x_{i}}{\partial\varphi_{j}} \\\frac{\partial y_{i}}{\partial T_{1}} & \frac{\partial y_{i}}{\partial T_{2}} & \frac{\partial y_{i}}{\partial T_{3}} & \frac{\partial y_{i}}{\partial\varphi_{j}}\end{bmatrix}} \\{= \begin{bmatrix}\frac{- f}{s_{x}Z_{i}} & 0 & \frac{- {fX}_{i}}{s_{x}Z_{i}^{2}} & \frac{- {f\left( {{Z_{i}\frac{\partial X_{i}}{\partial\varphi_{j}}} - {X_{i}\frac{\partial Z_{i}}{\partial\varphi_{j}}}} \right)}}{s_{x}Z_{i}^{2}} \\0 & \frac{- f}{s_{y}Z_{i}} & \frac{{- f}\; Y_{i}}{s_{y}Z_{i}^{2}} & \frac{- {f\left( {{Z_{i}\frac{\partial X_{i}}{\partial\varphi_{j}}} - {Y_{i}\frac{\partial Z_{i}}{\partial\varphi_{j}}}} \right)}}{s_{x}Z_{i}^{2}}\end{bmatrix}}\end{matrix}$

where f is the focal length of the C-arm fluoroscope; (s_(x), s_(y)) arethe respective x and y components of the C-arm fluoroscope pixel size;(X_(i), Y_(i), Z_(i)) are coordinates of a given point P on theestimated fiducial image generated in step 425, (T₁, T₂, T₃) are thetranslational elements, and φ is the rotational element.

Alternatively, the elements of the Jacobian Matrix corresponding topoints, lines, and ellipses may be respectively computed according tothe following relationships:

${\frac{\partial\left( {\delta \; x_{i}} \right)}{\partial P_{j}} = \frac{\partial x_{i}}{\partial P_{j}}},{\frac{\partial\left( {\delta \; y_{i}} \right)}{\partial P_{j}} = {\frac{\partial y_{i}}{\partial P_{j}}\mspace{14mu} ({point})}}$$\frac{\partial\delta_{i}}{\partial P_{j}} = {\frac{1}{\sqrt{A^{2} + B^{2}}}\left( {{A\frac{\partial x_{i}}{\partial P_{j}}} + {B\frac{\partial y_{i}}{\partial P_{j}}}} \right)\mspace{14mu} ({line})}$$\frac{\partial\delta_{i}}{\partial P_{j}} = {K\left( {{{\left( {{2{Ax}_{i}} + {By}_{i} + D} \right)\frac{\partial x_{i}}{\partial P_{j}}} + {\left( {{2{Cy}_{i}} + {Bx}_{i} + E} \right)\frac{\partial y_{i}}{\partial P_{j}}\mspace{14mu} ({ellipse}){where}\mspace{14mu} K\mspace{14mu} {is}\mspace{14mu} a\mspace{14mu} {scale}\mspace{14mu} {factor}\mspace{14mu} K}} = {\frac{{ab}\left( {{A\; C} - {B^{2}/4}} \right)}{\left( {2\sqrt{a^{2} + b^{2}}} \right)\det \; (M)}.}} \right.}$

The software may be designed such that the Jacobian Matrix is generatedany time between the determination of the current estimate and thebeginning of step 435. For example, in a multiprocessor data system 140,the Jacobian matrix and the Euclidean distance may be calculatedconcurrently on separate computers or processors. It will be apparent toone skilled in the art that different computing architectures may beemployed to minimize the time to perform steps 425-440.

In step 440, the software computes a new pose estimation based on thecurrent pose estimation, the Jacobian Matrix, and the Euclideandistance. Newton's method may be used to calculate the pose estimation,although other numerical optimization methods may be used, provided thatthe method is sufficiently robust and stable to compute a poseestimation in real time without suffering from singularities. In usingNewton's method, the software implements the following equation:

${{\Delta \; P} = \frac{\delta}{J}},$

where J is the Jacobian Matrix computed in step 430, δ is the Euclideandistance matrix computed in step 435, and ΔP is the change intranslation and rotation (six-DOF) between the acquired image and theestimated image. In a particular embodiment, the software computes andstores values for ΔP, and then segregates ΔP into a translationalcomponent and a rotational component. With these values, the softwarecomputes the new pose estimation as T_(i+1)=T_(i)+ΔT, (translation), andR_(i+1)=R_(i)ΔR_(i) (rotation), in which T_(i+1) is the translationalcomponent of the new pose estimation, T_(i) is the translationalcomponent of the current pose estimation, ΔT, is the translationalcomponent of ΔP, R_(i+1) is the rotational component of the new poseestimation, R_(i) is the rotational component of the current poseestimation, and ΔR, is the rotational component of ΔP. It may benecessary to normalize the Jacobian Matrix so that each degree offreedom has the same weight, and to normalize the Euclidean distancematrix so that each distance metric will have the same weight.Normalizing the Jacobian Matrix includes deriving a normalization matrixH that scales the parameter space without changing the least-squaressolution. The matrix H may be derived through the relationship

G(JH)(H ⁻¹ ΔP)=Gδ

where G is the Euclidean distance δx_(i), δy_(i), J is the JacobianMatrix, ΔP is the current pose estimation, and δ is the standarddeviation. The matrix H can be constructed by using a parametercovariance matrix, which requires an initial estimate of the solutionand the standard deviation of each parameter. The standard deviationsare derived by computing statistics across each iteration of steps425-445. For example, column scaling may be done using the mostconservative estimates for the H matrix elements. Column scaling is amethod known to the art, which does not use any a priori information,but requires identifiable parameters.

The matrix for normalizing the Euclidean distance matrix, G, may becomputed by using a Gauss-Markov estimation, which requires a goodestimate of the covariance matrix corresponding to the Euclideandistance. Another statistical approximation is to use a diagonal matrixwith each entry being the reciprocal of the standard deviation of aparticular element within the Euclidean distance matrix, which is oneapproach to make the Euclidean distance matrix a reasonable measure ofthe size of the least square error vector.

In step 440, the software computes the pose estimation using thenormalization matrices H and G according to the following relationship

${\Delta \; P} = {H \cdot {\left( \frac{G \cdot \delta}{G \cdot J \cdot H} \right).}}$

The software then stores the values corresponding to the new poseestimation.

In step 445, the software compares the Euclidean distance with apre-determined error threshold. If the Euclidean distance, whichcorresponds to the error in the pose estimation, is sufficiently small,the software stores the new pose estimation values as the final poseestimation. Otherwise, the software repeats the instructions associatedwith steps 425-440, with the new pose estimation taking the place of thecurrent pose estimation.

The “basic” algorithm 400 illustrated in FIG. 4 is exemplary, andvariations are possible. For example, once the Euclidean distance iscomputed in step 435, it may be compared with the pre-determined errorthreshold immediately. In this variation, if the Euclidean distance isbelow the pre-determined error threshold, the algorithm proceedsdirectly to step 450 without computing a new Jacobian Matrix and poseestimation.

FIG. 5 illustrates a “robust” pose estimation algorithm 500, which maybe implemented in step 315 as an alternative to the “basic” algorithmillustrated 400 in FIG. 4. The “robust” algorithm 500 is a two stageiterative process, wherein a plurality of course pose estimations arecomputed in a first iterative sub-process 560, and the pose with theminimum error is chosen to be the input pose estimation to a seconditerative sub-process 570. In the second iterative sub-process 570, aplurality of fine pose estimations are calculated, and the fine poseestimation with the least error is selected. The remaining plurality offine pose estimations may be used to provide a quality termcorresponding to variability of the final fine pose estimation.

In the “robust” pose algorithm 500, the software segments features 106on the fiducial device 105 from the acquired image in step 410 in amanner substantially similar to step 410 of the basic algorithm 400.

In the first iterative sub-process 560, in step 515 the softwareperturbs the current coarse pose estimation by some large perturbationfactor. If this is the first pass through the first iterativesub-process 560, the current course pose estimation may be an initialguess. Then the software computes a coarse pose estimation in step 520.The instructions executed in step 520 are substantially similar to thoseexecuted in the “basic” algorithm 400, except that only points are usedfrom among features 106 in computing the coarse pose estimation (asopposed to using lines and ellipses). In step 525, the software storesthe current coarse pose estimation values for later use and returns tostep 515, wherein the software perturbs the current coarse poseestimation, and repeats the estimation process of step 520. Theperturbation should be sufficiently large to provide a new currentcoarse pose estimation that is sufficiently distinct from the previouscoarse pose estimation such that the subsequent coarse pose estimationmay correspond to a distinct local minimum relative to the first coarseestimate. In an exemplary embodiment, a large perturbation involves atranslation of approximately 30 cm and a rotation of any angle. In otherwords, the “robust” algorithm 500 includes iterative coarse estimates toassure that the resulting set of coarse pose estimates encompass all ofthe possible local minima in the pose estimation space. The number ofiterations required for the first iterative sub-process 560 may beprovided by the operator, or may be a stored parameter value withinmemory 142 of data system 140.

In step 530, processor 141 retrieves the values for each of the coarsepose estimations that were stored in step 525, and selects the coarsepose estimation that corresponds to the smallest Euclidean distanceδx_(i), δy_(i). The selected coarse pose estimation is then used as aninput to the second iterative sub-process 570.

The software now performs second iterative sub-process 570. In step 540,the software perturbs the selected coarse pose estimation by a smallamount, and then uses the perturbed coarse pose estimation as theinitial pose estimation for step 545. For example, a small perturbationmay involve a translation of approximately 10 cm and a rotation ofapproximately 90 degrees. Step 545 may be substantially similar to steps425-450 of the “basic” algorithm 400, except that the initial poseestimate is the output of step 530. It is also similar to the softwareimplemented in step 520, except that step 545 may use some or all of thesegmented features 106 of the fiducial device 105, such as the points,lines, and ellipses. The result of step 545 is a current fine poseestimation, the values of which are stored for each iteration of thesecond iterative process 570. The software may repeat step 540-545,where each current fine pose estimation is perturbed in step 540 tobecome the input to step 545.

The number of iterations of the second iterative sub-process 570 may beinput by an operator through the user interface 145, or it may stored inmemory 142 as a configuration parameter value.

In step 550, processor 141 retrieves the stored fine pose estimationsand selects one to be the final pose estimation. For example, thesoftware may select the fine pose estimate corresponding to the minimumEuclidean distance δx_(i), δy_(i), and calculate the average of theEuclidean distances computed for each iteration of the second iterativesub-process 570. The computed average corresponds to a quality figure,which may provide the operator with feedback regarding the precision ofthe registration of the fluoroscopy image from first imager 115 to theTRUS image of second imager 110. Further, the software may “removeoutliers” whereby Euclidean distances that lie beyond a pre-defineduncertainty threshold are not included in the average. Also, in anotherexemplary embodiment, the software selectively retrieves a plurality offine pose estimations and their corresponding Euclidean distances, andcalculates a weighted average of the selected fine pose estimations,wherein the respective Euclidean distance corresponds to a weightingfactor or coefficient for the weighted average. It will be apparent toone skilled in the art that many algorithms for selection, outlierremoval, interpolation, and/or estimation may be implemented forselecting the fine pose estimation, all of which are within the scope ofthe invention.

Returning to exemplary process 300, regardless of whether the “basic”pose algorithm 400 or the “robust” pose algorithm 500 is used in step315, the result of step 315 includes stored values corresponding to apose estimation, which may be used to transform the image of thefeatures 106 on the fiducial device 105 according to six-DOF from thecoordinate frame of first imager 115 to the coordinate frame of fiducialdevice 105.

In step 320, the software registers the image from the first imager 115to the image of the second imager 110. If fiducial device 105 is fixedin position and orientation relative to second imager 110, then this maybe accomplished by the software, which transforms the entire (C-arm)image acquired by the first imager 115 into the coordinate frame ofsecond imager 110. However, there may be an intermediate transformationmatrix between the coordinate frame of fiducial device 105 to that ofsecond imager 110. If fiducial device 105 and second imager 110 arefixed relative to each other, the intermediate transformation matrix maybe constant, the values for which may be stored in the memory 142 ofdata system 140 as calibration data. However, second imager 110 may moverelative to fiducial device 105. In this case, the intermediate six-DOFtransformation matrix may change dynamically, and its values (positionand orientation) will have to be measured. The position and orientationof second imager 110 (the TRUS imager) relative to fiducial device 105may be measured by various encoding mechanisms such as optical encoders.Data system 140 may periodically retrieve these position and orientationmeasurements and store the corresponding values in memory 142. Oneskilled in the art will readily recognize how to implement anintermediate transformation matrix in the software.

Accordingly, the result of step 320 includes the image from first imager115 (C-arm) and an image from second imager 110 (TRUS) registered in asingle coordinate frame corresponding to the three-dimensional field ofview of second imager 110. The software then displays the C-armfluoroscopy image from first imager 115 superimposed over the TRUS imageinform second imager 110 in accordance with step 325, showing what isvisible to both imaging modalities within the field of view of secondimager 110.

Exemplary process 300 may repeat steps 305-330 in order to register anadditional C-arm fluoroscope image from first imager 115 to thecoordinate frame of second (TRUS) imager 110. In this case, it would bepreferable to have the next C-arm image from first imager 115 be from adifferent orientation relative to the previous C-arm image, so thatfeatures such as foreign objects 125 may be more easily located in thethree-dimensional TRUS coordinate frame of second imager 110. Forexample, in a first iteration of step 305, first imager 115 (C-arm) maybe reoriented into a position that may advantageously provide images ofthe foreign objects 125 in the surrounding soft tissue 120. It may bepreferable to reorient first imager 115 such that the subsequent firstimager 115 image coordinate frame is substantially orthogonal to theprevious first imager 115 coordinate frame. A preferable juxtapositionof first imager 115 for subsequent image acquisitions may depend on thenature of the subject and the medical procedure being performed. If theC-arm, or any imager being used for first imager 115, has multipledetectors and multiple fields of view, it may be possible tosimultaneously acquire multiple images from first imager 115 at multipleimage coordinate frames. Exemplary process 300 may be compatible withsimultaneous multiangle imaging. In this case, first imager 115 wouldonly need to be oriented once.

Subsequent iteration of steps 305-330 may be performed in a mannersubstantially similar to the above description for those steps. At thecompletion of step 330, the desired number of C-arm images from firstimager 115 have been acquired, registered according to their own uniqueorientation, and stored. Again, FIG. 2 illustrates two images 210 offeatures 106 on fiducial device 105 acquired by first imager 115 atmultiple angles. As illustrated, each image is at a differentorientation, and each provides a unique image of features 106.

FIGS. 6A and 6B show two different C-arm fluoroscopy images of a fieldof seeds 125 acquired at different viewing angles. The images shown inFIGS. 6A and 6B are exemplary of two images acquired and stored inexemplary process 200, except that the features 106 of the fiducialdevice 105 are not shown.

FIG. 7 illustrates a field of four seeds 125, all of which are withinthe fields of view of multiple orientations of first imager 115. Thecombination of unique orientations of the field of seeds 125, asprojected within the fields of view of first imager 115 at differentviewing angles, enables the use of one of several algorithms to matchand locate the seeds.

Returning to process 300, seed matching is performed in step 340.

FIG. 8 illustrates an exemplary process 800 for matching seeds in step340. As mentioned earlier, each seed is a small radiation implant thatis injected into a prostate for medical procedures such as cancertreatment. By knowing the location of each seed in relation to theprostate, it is possible to more accurately deliver a radiation dose tothe prostate, and more precisely estimate the radiation dose delivered.

In step 810, the software performs one of a number of known imageprocessing algorithms to identify and locate the seeds in each C-armfluoroscopy image from first imager 115. This results in a plurality ofsets of coordinates in the TRUS coordinate frame of second imager 110:one set per C-arm image from first imager 115, and one coordinate perseed. For a given set of seeds, each seed having been extracted fromeach C-arm fluoroscopy image, the software assigns a network flow. Inassigning a network flow, the software assigns a node to each seedidentified in each C-arm image and assigns a unique plurality of linksbetween each node in one of the C-arm images with a plurality of nodesin the other C-arm images. The plurality of links between the nodes ofeach C-arm image represents the possible matching of seeds between theC-arm images. The network flow may link nodes of multiple C-arm images.Nodes associated with three C-arm images may be networked accordingly.The result of step 815 includes a set of links covering all possiblecombinations of seeds.

In step 820, processor 141 executes the software to compute the “costmetric” of each link in the network. The cost metric refers to thelength of a link between nodes. Each possible matching of seeds has anaggregate cost metric, wherein the aggregate corresponds to the sum ofthe links between the nodes. The correct matching of seeds is that whichhas the least aggregate cost metric. The software may compute a costmetric for a given link by computing the least square 3D position of thenode such that the projected image coordinate of this node exhibits theleast deviation from its observed locations in the various images. Itwill be apparent to one skilled in the art that other approaches tocomputing the cost of a link are possible and within the scope of theinvention. For example, existing algorithms for computing minimum-costnetwork flow include Primal-dual algorithms for solving minimum-costflow problems.

In step 825, the software stores the matching of seeds with the minimalcost metric in memory 142. The stored values include an identifier for agiven seed and the 3-D coordinate for that seed. The result of thisapproach is a set of matched seeds that are identified across multipleC-arm fluoroscopy images from first imager 115. Seed matching, asperformed in steps 810-825, may be done in real time or near real time.

An alternative to the network flow of process 800 is to use theHungarian algorithm, which employs an assignment-problem-based approach,also known as the minimum weight bipartite matching problem. TheHungarian algorithm addresses the minimum-cost flow problem for twoC-arm images using a NxN cost matrix of identified seeds as input. Thedimension N pertains to the number of seeds identified in the two C-armimages. Implementations of the Hungarian algorithm are known to the art.

Exemplary process 800 described above may also be used to identify aparticular seed shape in addition to its location. Many seeds used inmedical procedures addressed by the present invention have varyinggeometrical features. The software may use the network flowimplementation described above to determine the shape of the seed byincorporating known information regarding the seed's geometry. Thesoftware may use the geometry in computing the cost metric by assigninga node to each of the endpoints as using the known seed parameters inevaluating the cost metric associated with the seed. For example, if theseeds are cylindrical in shape, the software may incorporate theendpoint locations information along with information regarding thelength of the seed. Accordingly, the cost computation function will useboth the end points of each seed and its length in three dimensions toevaluate the cost metric. Thus the process 800 may be fine-tuned toidentify seeds of a specific shape.

Variations to the above exemplary embodiments are possible and withinthe scope of the present invention. For example, fiducial device 105 maybe moved within the six-DOF space of second imager 110 while the firstimager 115 is left static. In this case, the movement of fiducial device105 may be tracked by using imagery from first imager 115. Bydetermining the locations of features 106 of the fiducial device 105 insix-DOF space relative to an anatomical feature of interest, andrepeating this for multiple positions of fiducial device 105, it may bepossible to determine the relative distances and orientations of(anatomical) objects in the six-DOF space. Surgical tools may be trackedusing a similar technique. This may also be useful in applications suchas tracking surgical tools.

In another variation, fiducial device 105 may be mounted on a roboticdevice (or any mobile/encoded mechanical actuator), in which theposition and orientation of fiducial device 105 in six-DOF space can beregistered to the coordinate frame of first imager 115, which in turnmay be registered to the coordinate frame of second imager 110.

Fiducial device 105 may also serve the additional purpose of onlinecalibration for first imager 115, wherein the need for preoperativecalibration can be alleviated either partially or completely using oneor more images of fiducial device 105 from either a single or multiplepose of first imager 115. Since the geometry of features 106 on fiducialdevice 105 are known, any distortions in the image of features 106provided by first imager 115 may be artifacts of spatial distortioncreated by first imager 115. By characterizing and compensating forspatial distortion in the images of the features 106, it is possible toprovide operating room personnel with imagery from first imager 115 inwhich spatial distortion is substantially mitigated. The products ofthis calibration procedure may include both the imager's spatial imagingparameters, which are generally constant, as well as variable distortioneffects, which are dynamic. Accordingly, fiducial device 105 may be usedfor spatial calibration in a preoperative manner, in a real time manner,or both.

Although the exemplary embodiment of the present invention involves anapplication in radiation treatment of a prostate, other applications mayinclude angiographic reconstruction of catheters and/or blood vessels.It will be readily apparent to one skilled in the art that many suchapplications are possible and within the scope of the invention.

In a variation to exemplary processes 400 and 500, steps 410(segmentation) and 425 (current pose estimation) may be iterated withinthe respective overall processes. Iteratively performing segmentationand pose estimation may not only improve accuracy but may also reducethe required accuracy of the segmentation step. This is becausesegmentation and pose estimation have an inherent duality, i.e. knowingone directly relates to the other. Accordingly, iteratively runningsteps 410 and 425 may improve the accuracy of the segmentation and thecurrent pose estimation with each pass. This may reduce or eliminate theneed to have an accurate initial segmentation.

FIG. 9 shows a second exemplary system 900 for registering ultrasoundimages and C-arm fluoroscope images according to the present invention.System 900 may include commercially available TRUS and C-arm fluoroscopyhardware. System 900 includes a C-arm fluoroscope 910; a C-armfluoroscope signal processor 915; a C-arm position controller 950; atransrectal sheath 920; a fiducial 925 disposed on the sheath; a TRUSprobe 930; a TRUS stepper 940; a TRUS probe position encoder 935; a TRUSsignal processor 945; a data system 955; and a computer 960, whichstores and executes software for implementing processes associated withthe present invention on data system 955. The computer 960 may includemultiple computers, including remote databases and embedded processors.It will be apparent to one skilled in the art that data system 955 andcomputer 960, may be provided in many different configurations.

FIG. 10 shows an exemplary process 1000 for registering ultrasoundimages and C-arm fluoroscope images according to the present invention.The process 1000 includes a C-arm fluoroscope pre-operation calibrationstep 1010; an image processing step (or image processing pipeline) 1020;a registration step 1030; a seed matching step 1040; and a dosimetrycomputation step 1050. Relevant aspects of the process 1000 areconsistent with the use of C-arm fluoroscope 910 in current clinicalprotocols.

During pre-operative calibration step 1010, calibration fixtures (X-rayfiducials) are attached to C-arm fluoroscope 910, and images acquired byC-arm fluoroscope 910 are processed by software running on data system955. The software then computes intrinsic spatial parameters of C-armfluoroscope 955.

Intra-operatively, system 900 uses C-arm fluoroscope 910 to acquireimages of a prostate. The software then processes the images in a mannersimilar to the processes used to compute the intrinsic spatialparameters. However, the software uses the intrinsic spatial parametersgenerated in step 1010 to produce spatially calibrated imagery. In doingso, the software may implement the same or similar processes andoperations as those of image processing pipeline step 1020.

Calibration of C-arm fluoroscope 910 in step 1010 may be performedpre-operatively, before the patient is brought into the room.Commercially available C-arms fluoroscopes generally apply an electronicimage intensifier, in which “parameter shift” may occur unexpectedly.Therefore it is recommend that calibration step 1010 be performed beforeeach procedure, even if the same C-arm fluoroscope 910 is being used. Itis assumed that intrinsic spatial parameters computed in step 1010 donot change during a single procedure.

Referring to FIGS. 12 a and 12 b during calibration step 1010, adewarping grid 1200 is placed on the detector of C-arm fluoroscope 910and C-arm signal processor 915 processes the images in a set of poses ina transverse and sagittal sweep planes at approximately +30°, 0°, −30°poses. Next, the calibration step 1010 uses an image processing pipelinesubstantially similar to that described for step 1020, only withoutusing the dewarping that is computed during calibration step 1010. Instep 1010, the software running on data system 955 assigns labels to thebeads of the dewarping grid 1200 automatically and derives a dewarpingfilter function. Next, a calibration object may be placed in the field,and the C-arm fluoroscope 910 is driven through the same poses asbefore. This time, however, the software running on data system software955 processes the images through the complete pipeline like that ofimage processing pipeline step 1020, now including dewarping. Thesoftware labels the calibration fixture automatically and thencalculates the intrinsic spatial parameters of C-arm fluoroscope usingtechniques known in the related art.

FIG. 12 a shows a dewarping grid structure 1200 that may universally fitC-fluoroscopes of different detector sizes. An adjustable band clamp 910goes around the C-arm fluoroscope 910, so that the plate 1220 hooks intoblocks positioned around the C-arm, as shown in FIG. 12 b. The plate1220 is firmly attached to image intensifier and substantially does notmove during the calibration process 1010. The dewarping process of imageprocessing pipeline step 1020 does not require that the centers of thedewarping plate 920 and the intensifier be aligned. The radio-opaquebeads 1230 may be arranged in triangular pattern. The dewarpingalgorithm implemented by the software running on data system 955includes three basic steps: (1) establish a correlation between pointsin the image and points in the known bead pattern 1230 by using a fastcrawling method that can circumvent problems caused by missing beads1230 in the image, (2) determine the rigid body transformation and scalefactor that will align the images, by using a conventional singularvalue decomposition method; and (3) align the images and performleast-squares polynomial fit between the point sets. The software mayuse nth-order Berstein polynomials to describe the patterns in thedewarping grid and then achieve deformable matching between warped andwarp-free features.

In step 1020, the software running on data system 955 acquires imagesfrom C-arm fluoroscope 910. The software then delineates and identifiesobjects (seeds, needles, fiducials, etc.) in the fluoroscopic images. Indoing so, image processing step 1020, includes a chain of operationsacting on the X-ray images acquired by C-arm fluoroscope 910: imageacquisition; dewarping; noise reduction; image normalization; backgroundremoval; thresholding; and pre-labeling.

In acquiring images from C-arm fluoroscope 910, C-arm signal processor915 may use hardware integral to fluoroscope C-arm 910. C-arm signalprocessor 915 may output digital signals or produce a VHS signal thatcan be captured through a frame grabbing card on the computer 960. Anexemplary frame grabber card includes a Matrox frame grabber, althoughothers may be used.

The dewarping operation helps remove the spatial distortion introducedby the image intensifier of the C-arm signal processor 915. The presentinvention includes a dewarping fixture and corresponding computerinstructions in the data system software 955, which are described later.

In the noise reduction operation, X-ray images acquired by C-armfluoroscope 910 may be averaged over multiple frames taken with C-armfluoroscope 910 in a single position. While this method substantiallypreserves image features, it also may increase radiation dose from C-armfluoroscope 910. For example, averaging five short-burst fluoroscopyimage frames may give excellent noise suppression, while exposing thepatient to a dose of one radiograph.

The image normalization operation includes the act of normalizing theimage for different attenuation rates corresponding to different objectsin an image. The intensity of X-ray energy falls exponentially as theray passes through the body, and the gray level of objects depends veryheavily on the thickness of other objects that the ray had to penetrateon its way. Logarithmic normalization of the image may result in simplesum of the respective attenuations of the seeds and body. As theimplanted seeds, needles, and fiducials tend to have large attenuationcompared to bone and soft tissues, they may stand out brightly inlogarithmic image.

The background removal operation is important to the accuracy of theexemplary process 1000 and its ability to find the centers andcenterlines of objects in the X-ray imagery. Calculating a simpleweighted centroid may not suffice, because non-point like regions mayskew the determination of the center and centerline along the gray levelgradient. An exemplary approach to the background removal operation mayinclude applying a linear gradient correction followed by Gaussianfiltering over a logarithmic normalized image. This approach may be usedfor implanted spherical markers and further expand the method for lineobjects and seeds. The background removal operation may use amathematical morphology approach known in the related art forfluoroscopic imaging of prostate implants by defining a chain of fourbasic morphological operations: erosion, dilatation, opening, andclosing.

The thresholding operation may not be done statically, because theintensity of seeds and objects also depends on the thickness of bonestructures. However, the optimal threshold may be found automatically byapplying entropy maximization over a two-dimensional histogram thatreflects both the gray levels and spatial relationship of pixels.

The pre-labeling operation is performed to delineate and identify allobjects other than seeds, then digitally remove those from the images.For example, the software running on data system 955 may use prior shapeinformation about the objects in the field of view. Two sorts of denseobjects may be expected in the field of view: implanted markers(needles, wire markers, spherical markers), and brachytherapy seeds.Needles and wire markers typically appear as thin lines that aretypically well segmentable by a Hough transform implementation. At anygiven time, only a very limited number of marker objects are generallypresent in an X-ray field. Knowing the shape, orientation, and estimatedrelative positions of these objects, the data system software 955 maylabel them automatically, or with minimal input from the operator. Ifthey are labeled once, the data system software 955 may remember wherethey are supposed to show up the next time and how they actually look inthe images. After a marker object is positively identified, the datasystem software 955 calculates its 3D locations for registrationpurposes and then digitally removes its silhouettes from the X-rayimages. This process may be repeated until only seeds are left in theimages that then undergo a pattern recognition process (seed matching)described below. It should be noted that some seeds may have beenobscured by marker objects and thus inadvertently removed from theimages.

The result of the above operations within image processing pipeline step1020 includes a distortion-free binary image, in which only artificialimplanted objects are present, while all details of the natural body(bone, ligaments, etc.) are removed from the picture. Subsequentlabeling can be further assisted by optimizing material density of thecustom-made objects (de-warp grid, calibration phantom, rectal sheath,fiducials, etc.) so that the processing pipeline step 1020 can befine-tuned to suppress or enhance these objects selectively.

In registration step 1030, the software running on data system 955cinoytes a sux-DOF transformation between the coordinate frame of C-armfluoroscope 910 and the coordinate frame of TRUS probe 930. In step1030, fiducial 925 is used as a reference that can be seen by both C-armfluoroscope 910 and TRUSprobe 930 and is substantially stationary withrespect to the prostate during process 1000.

Referring to FIGS. 13 and 14, the fiducial 925 may be coupled to thecoordinate frame of the TRUS probe 930 in many ways. An exemplary methodcan use a sheath 1320, either complete or partial, either inside oroutside the patient's body. In the event of using transcutaneousultrasound, fiducial 925 could be mounted elsewhere. Sheath 1320 couldbe used in addition to the mount to which fiducial 925 is attached. Anexemplary transrectal sheath 920 includes six-DOF external fiducialpattern 1310 mounted on the outside of the sheath 1320, that is mountedon the base of the TRUS stepper 940 of system 900. Inside sheath 1320the TRUS probe 930 can move freely without mechanical interaction withthe rectum, thus the location and shape of the prostate 1330 does notchange due to scanning motion of TRUS probe 930. Sheath 1320 may be madeof thin plastic, mylar, polypropylene, polyethylene or similar materialand may fit tightly round the TRUS probe 930.

As illustrated in FIG. 15A, gel material 1510 substantially providesultrasonic coupling between the sheath 1320 and the rectum wall. Thepresence of the sheath 1320 may cause some degradation of the acousticsignal, although minimizing the thickness of the sheath and maintainingeven distribution of the coupling gel 1510 mitigates this effect.

As illustrated in FIGS. 17A-17B, multiple air channels 1710substantially provide equalized air pressure in the sheath 1320 duringscanning motion of the TRUS probe 930.

Fiducial pattern 1310 on the outer surface of the sheath 1320 maycontain a suitable combination of straight lines, helixes, and ellipsesformed using an opaque substance, as illustrated in FIGS. 17A-17B. Inorder to determine the location of fiducials 1310 in the coordinateframe of C-arm fluoroscope 910, the relative pose of the C-arm imagesmust be known. Fiducial pattern 1310, which is preferablyprecision-machined, may be rigidly fixed within the field of view ofC-arm fluoroscope 910 to assist in the determination of the relativepose of the multiple images acquired by C-arm fluoroscope 910.

Fiducials pattern 1310 on sheath 1320 provides a common coordinatesystem for the individual X-ray images acquired by C-arm fluoroscope910. Further, sheath 1320 may be rigidly mounted on the base of TRUSstepper 940, and the motion of TRUS probe 930 may be encoded relative tothe base of TRUS stepper 940, so there is a known spatial relationshipbetween the TRUS probe 930 and fiducial pattern 1310 at any time.Therefore, fiducial pattern 1310 on sheath 1320 provides a commonreference coordinate system for both X-ray and the TRUS images. Thisframe of reference may be considered to be sufficiently stationary withrespect to the prostate during a registration session typicallyinvolving the acquisition of multiple X-ray shots and a TRUS volume.Fiducial pattern 1310 can be auto-segmented.

In seed matching step 1040, the software running on data system 955matches the radioactive seeds across multiple registered X-ray imagesacquired by C-arm fluoroscope 910. Reconstructing implanted seeds may becomplicated by a large number of seeds (sometimes as many as 150)confined in a small volume, overlapping one another in an apparentdisorder. The X-ray images may also be very noisy, especially in thelateral direction where the pelvic and femur bones create a strongheterogeneous shadow over the seeds. Accordingly, it is suggested thatthree X-ray images be acquired by C-arm fluoroscope 910. In addition,the apriori information on the inserted seeds can be used to help thetechnician position C-arm fluoroscope 910 to acquire images.

The software running on data system 955 may include instructions toimplement a technique that draws heavily on prior information fromTRUS-based implant plan (according to the surgeon's specifications) andfrom earlier seed matching and reconstruction sessions. Havingestablished 3D registration between the respective coordinate frames ofC-arm fluoroscope 910 and TRUS probe 930 in step 1030, the software maytransfer spatially registered information from the coordinate frame ofTRUS probe 930 to the coordinate frame of C-arm fluoroscope 910 toassist seed segmentation and matching in step 1040 information from theTRUS space to X-ray space to assist seed segmentation and matching, byrecognizing previously reconstructed seeds. The software may implementseveral methods for doing this.

Certain commercially available TRUS-based treatment planning systems,such as Interplant® made by CMS Burdette Medical Systems ofUrbana-Champaign, Ill, have the capability of spatially locating aninserted needle. Although there is some tissue deformation, a goodestimation is provided regarding where the seeds will end up in theprostate. By projecting the location of the needle forward onto the 2DX-ray pictures, a reasonably confined area may be defined where thesenew seeds are supposed to show up in X-ray. Consequently, the softwarerunning on data system 955 can label clusters of previously implantedold seeds that could not have come from the newly implanted needles. Thesoftware may update the 3D locations of old seeds and then excludesthose from seed matching step 1040. When a needle is retracted from theprostate, it leaves the seeds behind along a slightly bent trajectory.Based on observations of actual implants, the bending and perturbationof seeds from the ideal curvilinear trajectory may be statisticallymodeled by the software funning on data system 955. Using this estimatein confining the newly implanted seeds can improve the efficacy ofelimination of false candidates.

It may be assumed that seeds implanted by a previous fluoroscopy sessionhave not migrated substantially within the prostate since then. Byprojecting the location of these previously located seeds forward ontothe 2D X-ray pictures, the software running on data system 955 canobtain a reasonably confined area where they are supposed to show up inthe new X-ray images, then update their 3D locations and exclude thosefrom seed matching step 1040, as described above.

If the prostate did not move between imaging sessions, then bysubtracting the previous image from the new image taken by C-armfluoroscope 910 in a single pose may show the newly derived seedlocations and nothing else. In reality, the prostate moves somewhatbetween imaging sessions, and a pose of C-arm fluoroscope 910 isgenerally not precisely repeatable. A large portion of previouslyimplanted seeds can be excluded from the matching problem on this basis.The software can apply the images in processing pipeline step 1020 andperform deformable mapping using information that is mutual between theimages.

After seed matching step 1040 is complete, the software running on datasystem 955 may transform the location of seeds to the coordinate frameof TRUS probe 930, and their three-dimensional dose cloud may becomputed and analyzed by the Interplant® system. The number of seedsfound may be different from the actual number of implanted seeds, so theactivity of seeds may be adjusted so that the analysis shows thedistribution of truly implanted dose. Finally, the remainder of thetreatment plan is updated, in order to maintain to maintain optimal dosecoverage of the prostate. For visual evaluation of the registrationperformed in step 1030, fiducial pattern 1310 and other reconstructedstructures may be projected back to the coordinate frame of TRUS probe930 and superimposed in the 3D view provided by the Interplant® system.

As well as being used for intra-operative seed matching andregistration, process 1000 could be extended to post-operative care. Inthis approach, the last available implant plan could be projected ontopost-operative images and the individual seeds would be matched. Thesoftware running on data system 955 may repeat imaging sessionsregularly over an expanded period of time. As a result, it is possibleto track the path of each individual seed, which could provide aninsight to how the prostate changes shape, volume, and pose over time.Having collected data from a sufficiently large number of patients, thesoftware may predict statistical motion of seeds with respect to therelevant anatomy. This may prove to be helpful in understanding theprocess of post-operative edema, a crucially important issue in prostatebrachytherapy. The efficacy of projecting prior information to new X-rayimages can be greatly facilitated by synchronizing the insertionsequence with the poses of C-arm fluoroscope 910, so that the newlyinserted seeds are less likely obscure each other. This condition can beeasily met if the implant is executed by rows of needles and images areacquired by C-arm fluoroscope 910 after every row.

FIG. 18 illustrates an exemplary system 1800 for registering objects,such as brachytherapy seeds, in a 3-D image space using images acquiredat multiple poses of a single 2-D medical imaging system. System 1800may provide for tracking of a C-arm fluoroscope (or some other imagingmodality) and for reconstruction and registration of a plurality ofobjects, such as brachytherapy seeds, without the use of an externaltracker. System 1800 includes a C-arm fluoroscope 1805 having atransmitter 1810 and a detector array 1815, which are mounted to agimbal 1820. Detector array 1815 is connected to a processor 1825 thathas a memory 1830 encoded with software (hereinafter “seedreconstruction software”) that implements processes associated with thepresent invention. User interface 1835 is connected to processor 1825.

C-arm fluoroscope 1805 acquires images of an anatomical region ofinterest 1845 within a patient 1840. Implanted within anatomical regionof interest 1845 is a plurality of seeds 1850. Using gimbal 1820, C-armfluoroscope acquires multiple fluoroscopy images of anatomical region ofinterest 1845 at multiple orientations or poses.

FIG. 19 illustrates three exemplary images acquired of seeds 1850acquired at a first pose 1905 a, a second pose 1905 b, and a third pose1905 c. Gimbal 1820 orients transmitter 1810 and detector array 1815 toprovide an image of seeds 1850 at each pose 905 a-c. As illustrated inFIG. 19, pose 1905 b differs from pose 1905 a by a rotation R₁ and atranslation T₁; and pose 1905 c differs from pose 1905 b by rotation R₂and translation T₂.

FIG. 20 is a diagram of an exemplary process 2000 for registeringobjects in a 3-D image space using images acquired at multiple poses ofa single 2-D medical imaging system. Process 2000 may be performed allor in part by the seed reconstruction software running on processor 1825in conjunction with user interaction via interface 1835. At theinitiation of process 2000, C-arm fluoroscope ready to acquire images ofanatomical region of interest 1845 (and seeds 1850) within the field ofview of detector array 1815. Further, an initial pose estimation ofC-arm fluoroscope 1805 may be provided and stored in memory 1830. Theinitial pose estimation may be an estimated rotation and translationfrom some default baseline orientation and position.

In step 2005, C-arm fluoroscope 1815 acquires an image of anatomicalregion of interest 1845 containing seeds 1850. In doing so, the user mayinput a command via user interface 1835. The seed reconstructionsoftware receives the command and in turn issues the appropriateinstructions to transmitter 1810 and detector array 1815 to acquire animage. The seed reconstruction software receives image data fromdetector array 1815 and stores the data values in memory 1830.

In step 2010, the seed reconstruction software segments the imageacquired in step 2005. In doing so, the seed reconstruction software mayimplement any of a number of segmentation image processing algorithmsthat are capable of identifying the pixels that correspond to seeds1850. The seed reconstruction software stores the segmented seedinformation in memory 1830.

In step 2015, the seed reconstruction software determines if anotherimage at another C-arm image (at another pose) is to be acquired. Forexample, three C-arm images, each at a different pose 1905 a-c may beacquired and used to reconstruct seeds 1850 and register them to a 3-Dimage space. More or fewer images may be acquired and used. Althoughmore C-arm images, acquired at different poses, may increase theaccuracy and precision of the reconstruction and registration of seeds1850, this occurs at the expense of computation time. In general, threeimages at poses 1905 a-c may strike a balance between accuracy/precisionand computation time. Variations in the number of images and poses mayvary, depending on the complexity and geometry of seeds 1850 and thecomputational power of processor 1825. One skilled in the art willreadily appreciate that such variations to process 2000 are possible andwithin the scope of the invention.

If another image is to be acquired at another pose, process 2000proceeds via the “yes” branch of step 2015 and proceeds to step 2020.

In step 2020, C-arm fluoroscope 1805 is positioned at a new pose. To dothis, the seed reconstruction software may issue commands to an actuator(not shown) coupled to gimbal 1820 of C-arm fluoroscope 1805 to achievea new pose. Alternatively, a user may command the actuator to a new posevia user interface 1835.

In step 2025, the new pose is estimated. If the seed reconstructionsoftware issued the command to C-arm fluoroscope 1805 in step 2020, thenew pose estimate may be the commanded pose. Similarly, if the userprovided the command to C-arm fluoroscope 1805 via user interface 1835,then the user-provided command may be the new pose estimate. In eithercase, the seed reconstruction software may store the estimated pose, andthe estimated transformation from the previous pose to the present pose.Referring to FIG. 19, the transformation from pose 1905 a to pose 1905b, for example, may be represented as rotation R₁ and translation T₁Rotation R₁ may be represented as rotation matrix or a quaternion.Translation T₁, may be represented as a vector. The pose estimation,which includes rotation R₁ and translation T₁ may be a rough estimate.Accordingly, C-arm fluoroscope 1805 does not require an externaltracker.

Having established a new pose for C-arm fluoroscope 1805 in step 2020,and having stored pose estimates in step 2025, the seed reconstructionsoftware may repeat steps 2005 and 2010. Accordingly, after multipleiterations of steps 2005-2025, multiple images are stored in memory1830, each of which having a corresponding pose estimation.

Returning to step 2015, if no more images are to be acquired, process2000 proceeds via the “no” branch of step 2015 and proceeds to step2020.

In step 2030, the seed reconstruction software computes a minimum costcorrespondence between the segmented seeds imaged in steps 2005 andsegmented in step 2010. The seed reconstruction software may do so byfinding the minimum of the following relation:

$\begin{matrix}{\arg \; \min} \\{R,T,M,f}\end{matrix}{\sum\limits_{i = 1}^{N_{1}}\; {\sum\limits_{j = 1}^{N_{2}}\; {\sum\limits_{k = 1}^{N_{3}}\; {C_{i,j,k}f_{i,j,k}}}}}$

where N₁, N₂, N₃ respectively refer to the number of seeds segmented inthe first, second, and third image; R is a combined rotation matrix,which includes rotations R₁ and R₂; T is a combined translation vector,which includes translations T₁ and T₂; C_(i,j,k) is a cost factor ofmatching the ith seed segmented in the first image with the jth seedsegmented in the second image with the kth seed segmented in the thirdimage; and M is the projection model of C-arm fluoroscope 1805. Theprojection model M refers to the spatial calibration parameters of C-armfluoroscope 1805. The projection model M is typically approximated as afive-parameter projection camera model, which may be assumed constantfor purposes of the present invention.

The term f_(i,j,k) is a binary variable that denotes a correspondencebetween the the ith seed segmented in the first image with the jth seedsegmented in the second image with the kth seed segmented in the thirdimage. In other words,fis a matrix that maps each seed segmented in oneimage with its counterpart in the other images.

In computing the minimum cost correspondence in step 2030, the seedreconstruction software iteratively selects a given correspondencef_(i,j,k), and then computes its associated cost factor C_(i,j,k). Theassociated cost factor C_(i,j,k) is a function of the combined rotationmatrix R and the combined translation vector T, both of which arederived from the rotations R₁ and R₂. and translations T₁ and T₂estimated in step 2025.

The seed reconstruction software may compute the cost factor C_(i,j,k)associated with a given correspondence f_(i,j,k) by computing aprojection error for the given correspondence. The seed reconstructionsoftware may do this by generating a line from the segmented seed ineach image to transmitter 1810 of C-arm fluoroscope according toprojection model M. The seed reconstruction software then computes anintersection point in 3-D space at the intersection of these generatedlines. Then, the seed reconstruction software projects this 3-D pointonto detector array 1815 according to the pose estimations determined instep 2025. Then, for each projection, the seed reconstruction softwarecomputes the distance from the projected point in the correspondingimage to the segmented seed in that image. In doing so, the seedreconstruction software computes three error distances, one for eachpose and corresponding image, which it may average to compute costfactor C_(i,j,k).

The seed reconstruction software may compute the cost factor C_(i,j,k)for each segmented seed (the ith seed segmented in the first image withthe jth seed segmented in the second image with the kth seed segmentedin the third image) for all possible correspondences f_(i,j,k). However,given all the possible combinations, this may be considered an “NP-Hard”combinatorial optimization problem. To improve computational efficiency,the seed reconstruction software may use a network flow-basedcombinatorial optimization solution, as discussed above.

Further to step 2030, the seed reconstruction software computes thesummed cost factors C_(i,j,k), resulting in a summed cost C for eachcorrespondence f_(i,j,k). The correspondence associated with the lowestsummed cost C is the minimum cost correspondence.

In step 2035, having solved for the minimum cost correspondencef_(i,j,k) in step 2030, the seed reconstruction software may store theminimum cost correspondence f_(i,j,k) in memory 1830.

In step 2040, the seed reconstruction software computes a minimum costpose estimation corresponding to the minimum cost correspondencef_(i,j,k) computed in step 2035. In doing so, the seed reconstructionsoftware may select a plurality of combined rotation matrices R and aplurality of combined translation vectors T, each of which have valuesthat vary from those of the pose estimation generated in step 2025. Thenumber of values, the range by which the values vary from the poseestimation generated in step 2025, and the spacing of the values, may bestored in memory 1830 as configuration parameters. These configurationparameters may be provided by the user via user interface 1835.

In step 2040, the seed reconstruction software computes a minimum costpose estimation, by which it computes the combined rotation matrix R andcombined translation vector T according to the relation described abovewith respect to step 2030

$\begin{matrix}{\arg \; \min} \\{R,T,M,f}\end{matrix}{\sum\limits_{i = 1}^{N_{1}}\; {\sum\limits_{j = 1}^{N_{2}}\; {\sum\limits_{k = 1}^{N_{3}}\; {C_{i,j,k}{f_{i,j,k}.}}}}}$

In this step, correspondence f_(i,j,k) is set to that computed in step2030, and the combined rotation matrix R and combined translation vectorT are varied. For each iteration of R and T, the seed reconstructionsoftware may compute the associated cost factor C_(i,j,k) of the ithseed segmented in the first image with the jth seed segmented in thesecond image with the kth seed segmented in the third image. The seedreconstruction software may do so by computing a projection error in thesame or similar manner as described with regard to step 2030 above. Foreach iteration of R and T, the seed reconstruction software computes asummed cost factor C according to the relation above. The seedreconstruction software then selects the pose estimation associated withthe lowest summed cost factor C. Step 2040 results in a minimum costpose estimation, which is in the form of combined rotation matrix R andcombined translation vector T, which in turn corresponds to a bestestimation for rotations R, and R2 and translations T1 and T₂,associated with the correspondence f_(i,j,k) computed in step 2030.

In step 2045, the seed reconstruction software determines if the poseestimation computed in step 2040 is optimal. In doing so, the seedreconstruction software compares the lowest summed cost C computed instep 2040 with an expected minimum cost. The expected minimum cost maybe stored in memory 1830 as a configuration parameter. Alternatively,the seed reconstruction software may display the lowest summed cost Ccomputed in step 2040 to the user via user interface 1835. In this case,the user may decide if the minimum cost pose estimation computed in step2040 is optimal, or if the process should be repeated with a newcorrespondence f_(i,j,k)

If it is determined that the pose estimation computed in step 2040 isnot optimal, then process 2000 proceeds via the “no” branch of step 2045to repeat steps 2030-2040. In the subsequent execution of step 2030, theseed reconstruction software computes a new minimum correspondence (step2030) corresponding to the new pose estimation computed in step 2050.With this new minimum correspondence, the seed reconstruction softwarecomputes a new minimum cost pose estimation (step 2040). Accordingly,the loop of steps 2030-2045 may be iterated whereby each iteration mayconverge on an optimal solution for the correspondence and the poseestimation.

If it is determined that the pose estimation computed in step 2040 isoptimal, then process 2000 proceeds via the “yes” branch of step 2045 tostep 2050.

In step 2050, the seed reconstruction software stores the minimum costpose estimation computed in step 2040 in memory 1830.

In step 2055, the seed reconstruction software reconstructs theplurality of seeds 1850 in a three-dimensional space using the first,second, and third segmented images, along with the minimum costcorrespondence and minimum cost post estimation between the images. Indoing so, the seed reconstruction software may assign x, y, and zcoordinates to each of the segmented seeds 1850 in a coordinate space.The coordinate space may be referenced to the optical axis of C-armfluoroscope 1805, although other coordinate references may be used.

In step 2060, the seed reconstruction software may register thereconstructed three-dimensional segmented seeds 1850 to athree-dimensional image space corresponding to an image acquired by athree-dimensional imaging modality, such as an ultrasound imager, andthe like.

Variations to system 1800 and process 2000 are possible. For example, ifthe plurality of seeds 1850 are insufficient to enable process 2000 toyield a sufficiently precise minimum cost pose estimation, thenadditional objects may be introduced, such as a wire or metal beads, byattaching them to the skin of patient 1840, and within the field of viewof detector array 1815 for each of poses 1905 a-c. Further, ifnecessary, a fiducial device may be introduced having a feature, such asan ellipse of fiducial device 105 illustrated in FIG. 1.

Further, although the descriptions of system 1800 and process 2000involve a C-arm fluoroscope, one skilled in the art will readilyappreciate that the reference to a C-arm fluoroscope is exemplary, andthat any two-dimensional medical imaging modality may be used providedthat the imaging modality enables image acquisition of anatomical regionof interest 1845 from multiple orientations or poses.

Additionally, although system 1800 and process 2000 involve imaging,segmenting, and reconstructing seeds 1850, such as brachytherapy seeds,one skilled in the art will readily recognize that other objects may besubstituted for the seeds 1850 described herein. Other objects mayinclude surgical tools, foreign objects, or any other plurality ofobjects within anatomical region of interest 1845 within patient 1840,wherein the objects are visible to C-arm fluoroscope 1805 (or otherimaging modality, as stated above) with sufficient contrast to besegmented using known segmentation techniques.

It will be apparent to those skilled in the art that variousmodifications and variation can be made in the present invention withoutdeparting from the spirit or scope of the invention. Thus, it isintended that the present invention cover the modifications andvariations of this invention provided they come within the scope of theappended claims and their equivalents.

1-10. (canceled)
 11. A device for registering images acquired by aplurality of medical imaging modalities comprising: a support structure;a first wire disposed on the support structure, the first wiresubstantially formed in the shape of a straight line; and a second wiredisposed on the support structure, the second wire substantially formedin the shape of an ellipse.
 12. The device of claim 11, wherein thesupport structure includes a hollow cylinder.
 13. The device of claim11, wherein the first wire and second wire include Tantalum.
 14. Thedevice of claim 11, wherein the first wire includes a diameter of about1 mm.
 15. The device of claim 11, wherein the second wire includes adiameter of about 1 mm. 16-27. (canceled)